Image data from the years 1959/1960 and 1999/2000 reveal a 2.4% decrease in the surface area of the Devon Ice Cap, Nunavut, over the last 40 yr. This has resulted primarily from extensive retreat of tidewater glacier margins on the eastern side of the ice cap, and shrinkage of its near-stagnant southwestern arm. Thinning of the ice cap has also increased bedrock exposure in the ice cap interior. However, since 1960 the northwestern margin of the ice cap has advanced slightly. Volume loss associated with these changes was estimated at −67 ± 12 km3 as calculated from two independent techniques. A digital elevation model (DEM) of the ice cap surface was used to delineate interior ice divides allowing patterns of change to be investigated at the drainage basin scale. Strong correlation between the hypsometric characteristics of drainage basins and the observed changes in ice-cap geometry suggests that these changes reflect interbasin differences in the inherent sensitivity of glacier mass balance to recent climate forcing. Response time calculations indicate that most of the ice cap is responding to recent climate warming, whereas the northwestern region is likely still responding to cooler conditions that prevailed during the Little Ice Age (LIA).
Introduction
Global climate models consistently predict that anthropogenic climate warming will be manifested most strongly in northern high latitudes (Mitchell, 1995; Johns et al., 1997; IPCC, 2001). There is already considerable evidence that predicted changes are occurring (Overpeck et al., 1997; Serreze et al., 2000). The predicted temperature changes would impact strongly on the mass balance and extent of the glaciers, ice caps, and ice sheets in the region, and the resulting shrinkage of these ice masses could contribute significantly to global sea-level change over the next century (Meier, 1984; Van de Wal and Wild, 2001; Meier and Dyurgerov, 2002).
Outside Greenland, the largest amounts of land ice (∼150,000 km2) in the Northern Hemisphere are found in the Canadian Arctic islands (Koerner, 1966). However, little is known about recent changes in the extent of these ice masses. Areal change measurements have been made for a section of the Barnes Ice Cap on Baffin Island (Jacobs et al., 1996) and for the small, stagnant Murray Ice Cap on northern Ellesmere Island (Braun et al., 2001). Over the period 1961 to 1993, the area of the measured sector of the Barnes Ice Cap decreased by ∼1%, while that of the Murray Ice Cap decreased by ∼28% between 1959 and 2000. Here we report results of an investigation into recent changes in the extent of the Devon Ice Cap, Nunavut. Investigation of glacier change at the scale of the whole ice cap is of particular interest because there are strong spatial gradients in climate, mass balance characteristics, glacier geometry, and terminus conditions across the ice cap. Thus, different sectors of the ice cap may be subjected to different climate forcings, they may have different inherent mass-balance sensitivity to specific forcings, and they may respond to these forcings at different rates.
Glacier response times are influenced by a combination of ice-mass geometry, mass-balance gradient, and terminus ablation rates (Johannesson, 1989; Bahr et al., 1998). All of these factors vary significantly across the Devon Ice Cap. Differences in mass-balance gradient between the east and west sides of the ice cap have been attributed to the effect of the Baffin Bay moisture source on regional accumulation patterns. These patterns are believed to have influenced the long-term evolution of the ice-cap geometry such that glaciers draining to its eastern margin occupy larger basins and descend to lower elevations than glaciers draining to the west (Koerner, 1977a). Asymmetry of mass-balance gradient and ice-cap geometry resulting in large differences in terminus ablation rates between the east and west margins may have a significant effect on response times characteristic of these regions.
This study quantifies recent changes in the areal extent of the Devon Ice Cap by comparing the position of ice margins extracted from 1959/1960 aerial photographs and 1999/2000 satellite imagery. Information on the topography of the ice-cap surface, obtained from the Canadian Digital Elevation Dataset (CDED), was used to delineate drainage basins within the ice cap, and changes were analyzed at the scale of individual drainage basins. Relationships among the magnitudes of observed changes, the characteristic response times, and the hypsometry of drainage basins were investigated to provide insight into the mechanisms controlling ice margin variations. Also, two independent methods were used to estimate volume change for selected basins and for the ice cap as a whole.
Study Area
The Devon Ice Cap occupies approximately 14,400 km2 (between 74°30′N and 75° 50′N and 80°00′W and 86°00′W) on Devon Island, which is located in the southeast of the Queen Elizabeth Islands, Nunavut, Canada (Fig. 1). Recent radio echo sounding data indicate that the current ice-cap volume is 4110 ± 140 km3 with a maximum ice thickness of approximately 880 m at the head of the eastward flowing basins (Dowdeswell et al., submitted). The highest elevation is 1921 m at the ice-cap summit. The eastern margin of the ice cap faces the North Open Water (NOW) polynya at the head of Baffin Bay. Large outlet glaciers that descend to sea level in this region experience high rates of accumulation (up to 500 kg m−2 a−1) and surface ablation (up to 2000 kg m−2 a−1) as well as mass loss by iceberg calving, resulting in relatively steep mass-balance–elevation gradients (2.7 kg m−2 a−1 m−1; Koerner, 1970). Smaller tidewater glaciers reach the ocean along the northern and southern margins of the ice cap. By contrast, the western margin of the ice cap terminates entirely on land at elevations of 300 to 500 m. This part of the ice cap is located in a precipitation shadow, where rates of accumulation and ablation (up to 200 kg m−2 a−1 and ∼−1000 kg m−2 a−1, respectively), along with mass-balance gradients (1.5 kg m−2 a−1 m−1) are much reduced. An arm extending approximately 80 km to the southwest of the main ice cap consists largely of near-stagnant, ablating ice. Nearly this entire sector of the ice cap lies below the equilibrium-line altitude (ELA) as defined by mass-balance data compiled by Koerner (1970).
Proxy data and field measurements provide insight into the climate history of the ice cap since the end of the Little Ice Age (LIA). Records of summer melt layers from ice cores in the Canadian Arctic indicate an abrupt warming trend beginning ca. AD 1850, with a second significant temperature increase after ca. 1925 (Koerner, 1977b). Stake measurements made across the northwest part of the ice cap (Sverdrup Glacier) since 1960 indicate a negative mass balance since observations began, with a trend towards increasingly negative balances since the late 1980s (Koerner and Lundgaard, 1995).
Methods
Changes in the surface area of the Devon Ice Cap were determined from remotely sensed imagery acquired in 1959/1960 and 1999/2000. Ice margins were digitized for each data set and drainage basin boundaries common to both data sets were delineated. Overlay analysis and raster cartographic techniques were then employed to detect and quantify areal differences between the extracted data sets.
DATA SOURCES
The satellite image data used in this study were obtained from 3 Landsat 7 ETM+ panchromatic (15-m resolution) scenes acquired on 13 July1999 (path 38 rows 6 and 7) and 2 August 2000 (path 36 row 7). All satellite imagery was purchased as L1G (radiometrically and geometrically corrected) processed data. The aerial photography consisted of 250 1:60,000 photographs acquired in late July and early August of 1959/1960 by the Government of Canada. The digital elevation model (DEM) used in this study was a subset of the Canadian Digital Elevation Dataset (CDED) produced from the National Topographic System (NTS) 1:250,000 map sheets which, in turn, were derived from the 1959/1960 aerial photography mentioned above. The DEMs were re-projected from the WGS-84 geographic coordinate system (3 × 6 arcsec) to a 100 m resolution NAD83 UTM grid. Vertical accuracy of the DEM is ±20 m over bedrock and ice margins decreasing to ±50 m throughout the interior regions of the ice cap (Gagne, pers. comm., 2002).
IMAGE PREPARATION
All image data were referenced to the UTM projection on the NAD83 datum. Satellite images were georeferenced to the 1:250,000 NTS map sheets using at least 40 ground control points (GCPs) that were clearly identified in both images (RMS error <60 m). The satellite imagery was orthorectified in PCI OrthoengineTM V8.0 using the CDED DEM to correct for terrain distortions. The orthorectified images were manually mosaiced to produce a 15-m resolution orthoimage of the entire Devon Ice Cap. The aerial photographs were obtained as contact prints and digitized at 300 dpi (5 m ground resolution) using a Vidar 36″ upright scanner. Each digital photograph was georeferenced to the Landsat 7 orthomosaic by selecting between 7 and 12 GCPs (depending on the amount of exposed bedrock) over bedrock features clearly identifiable in both images (RMS error <15 m). The digital photos were then orthorectified in the PCI OrthoengineTM V8.0 using the CDED DEM to correct for terrain distortions. Four mosaics of the aerial photography were generated in OrthovistaTM (image mosaicing software) and clipped to the geographic extents of the corresponding 1:250,000 map sheets.
LINEWORK CAPTURE
Ice margins and boundaries of interior bedrock regions were captured through on-screen digitizing techniques in ArcView 3.2TM geographic information system from both orthorectified digital data sets. Linework representing ice margins was digitized in “point mode” as a series of line segments and later converted to single polygon coverage in ArcInfo 8.0TM. Interior bedrock regions were digitized directly as polygons in ArcView 3.2TM. The 1960 ice margin was digitized from the contiguous ice mass. Any parts of the ice cap that were separated from the main ice cap by melting between 1960 and 1999 were included in the 1999 surface area calculations.
DRAINAGE BASIN DELINEATION
Drainage basins were delineated to allow areal changes to be associated with distinct catchment regions (Fig. 2). A river drainage network for the ice-cap surface was initially simulated from the DEM and overlaid on the satellite image mosaic.
Using the “Basin” utility in ArcView 3.2TM, outlet locations (points where the simulated drainage paths flow off the ice-cap edge) were identified interactively and the associated catchment areas were delineated based on flow accumulation and direction grids. Due to inaccurate representation of the ice surface by the DEM in some areas, particularly towards the margins, ice flow (as observed from the Landsat image) is not always consistent with the downslope direction derived from the DEM. In these cases, manual editing of the boundaries was necessary to prevent drainage basin boundaries from cross-cutting ice-flow lines. This was done through visual interpretation of the Landsat orthomosaic. Drainage basin boundaries were then clipped by both the 1960 and 1999 ice-margin polygons.
AREA CALCULATIONS
Calculation of area change in each drainage basin was performed using raster overlay techniques in ArcView GIS. The drainage basin boundary layer was clipped by the 1960 and 1999 ice-margin polygons to produce the same internal drainage divide structure for both years. The polygon coverages and interior bedrock regions were then converted to a 15-m resolution raster grid. Cells were coded with unique arbitrary values assigned to each polygon with the same polygons having the same coded value for both years. For each data set (1959/1960 and 1999/2000) the interior bedrock regions were subtracted from the ice surface polygon coverage, resulting in “ice surface only” coverages. Areal changes were then calculated by comparing tabulated pixel counts for each gridded ice-surface coverage on a per polygon basis.
ERROR ANALYSIS
Errors associated with the ice-cap area measured from the air photos were a function of both the accuracy of the coregistration between the air photos and the Landsat 7 ETM+ orthomosaic and the accuracy with which the ice-cap margin was identified and digitized. Errors of the Landsat 7 ETM+ derived areas were a function of the accuracy with which the ice-cap margin was identified and digitized only. Ice-cap margins were digitized using on-screen digitizing in ArcView GIS that allowed the operator to zoom in to the pixel level for point data collection. Digitizing inaccuracies over unobscured sections of the ice-cap margin were therefore based on (a) image resolution and (b) contrast between the ice-cap surface and the adjacent terrain.
Coregistration errors between the air photo orthomosaics and the Landsat 7 ETM+ orthomosaic were estimated to be ±22 m based on residual values between 40 (10 points from each air photo mosaic) independent checkpoints identified in each data set. However, for the area measurement to be affected by the maximum residual distance, the coregistration error would have to displace each mosaic in equal and opposite directions from the other mosaics. Since it is more likely that these errors are random, a value of ±11 m (50% of the coregistration error) was chosen as a more realistic estimate of the error attributed to coregistration (see Table 1).
Measurement uncertainty of the unobscured sections of the ice margin was judged to be no greater than the pixel size of the image being digitized (±15 m for the satellite mosaic and ±5 m for the aerial photography mosaic). This error was calculated as the product of the total ice-cap boundary length minus the length obscured, and the pixel width. For obscured sections of the ice-cap perimeter, error was estimated as the product of the obscured segment length and a maximum estimated line offset value applied to each type of obscurity (Table 1). Because all missing or obscured sections of the aerial photography orthomosaic were digitized from the 1:250,000 NTS map sheets, estimated offset distances were based on the coregistration accuracy between the aerial photography orthomosaic and the 1:250,000 NTS map sheets (±75 m).
Offset distances were estimated for the satellite orthoimage based on the confidence with which the ice margin could be identified within each type of obscurity. Because the margin was generally not completely obscured by cloud and shadow, error along these sections was estimated at ±90 m. Sections of the ice margin under cloud were faintly (intermittently) identified through semitransparent (broken) cloud cover. Similarly, sections of the ice-cap margin obscured by shadow were enhanced by applying a contrast stretching function to the imagery that allowed the margin to be identified with some degree of confidence in these areas. Ice margins obscured by late-lying snowpack, however, were less distinct and identification accuracy was estimated at ±120 m. Along these sections, the margin was digitized at the angle of inflection between the ice-cap boundary and the snow-covered ground that was highlighted by differential illumination as a result of the local incidence angle of the sun. In total, the ice-cap margin was obscured (either partially or fully) by 3% and 0.8% in the satellite imagery and the aerial photography, respectively.
Additional error could potentially have been generated in the orthorectification process due to registration errors between the imagery and the DEM. The DEM used in this study is essentially a digital version of the NTS topographic map sheets to which the Landsat 7 imagery was georeferenced with an RMS error of ±75 m. The aerial photography was subsequently georeferenced to the Landsat 7 imagery (RMS error of ±15 m) increasing the registration error between the aerial photography and the DEM to ±90 m. Therefore, the uncorrected aerial photography and Landsat 7 imagery may have been orthorectified based on topographic information offset by up to 90 m from the associated features on the imagery. Orthorectification errors resulting from inaccurate coregistration between the DEM and imagery however were not considered in this study because both data sets would be affected equally by this problem therefore the relative difference in area measurements would be minimal.
AREA ESTIMATES AND TOTAL ERROR
The total error (TE) of the area measured from each data source was computed from the formula;
where g1 … gn represent the individual sources of error as specified in Table 1. The total area (TA) and the associated error as measured from the 1999 Landsat imagery therefore is; and for the 1960 aerial photography is;
The error associated with the change of areal extent (difference between the 1999 area and the 1960 area measurement) was taken as the greater of the two error estimates.
Therefore,
The minimum error of the area change measurements derived for an individual drainage basin was estimated to be accurate to within ±11% for basins with unobscured margins. For all other drainage basins, error was estimated as a function of the length (and type) of the obscured margin within the basin.
The error analysis presented in this study likely overestimates the actual error because it assumes a maximum digitizing offset along the entire perimeter of the ice-cap boundary. Previous work based on repeat digitizing experiments has demonstrated that errors associated with digitizing points along the perimeter of a polygon are limited to a statistically derived epsilon band which is significantly smaller than the maximum offset width (Dunn et al., 1990). It is beyond the scope of this study however, to minimize the stated error by deriving a statistical distribution of digitized points about the ice-cap margin.
Results
The Devon Ice Cap has experienced a net decrease in area of 338 ± 40 km2 over the past 40 yr. Surface area changes were most significant in four distinct regions of the ice cap.
Retreat of the four major tidewater glaciers on the east coast (Fig. 3) accounted for 38 ± 4 km2 of the decrease in area of the ice cap. The most northerly glaciers (i) and (ii) experienced retreat of up to 1300 m from the 1960 terminus location and the central glacier (iii) retreated approximately 750 m. The most southerly glacier (iv) experienced maximum retreat of up to 3000 m.
Increase in the area of exposed bedrock in the interior regions of the ice cap accounted for 45 ± 5 km2 of the decrease in ice surface area (Fig. 3). The basins affected the greatest are located in the mountainous northeast and southeast regions of the ice cap. The increase in exposed bedrock area in these regions indicates a lowering of the ice-cap surface.
Marginal retreat of the southwest arm (Fig. 4) accounted for a 200 ± 17 km2 decrease in surface area of the ice cap.
Ice margin advance by an average of 130 m along a 80 km section of the northwest margin (Fig. 5) accounted for a 10.5 ± 2 km2 increase in surface area of the ice cap.
The changes described above account for approximately 76% of the total change experienced by the Devon Ice Cap. The remaining changes in terms of area decrease have occurred primarily as a result of glacier retreat within small basins located throughout the south and northeast sectors of the ice cap. The remaining increases in area have resulted from advance of the Sverdrup Glacier (basin 12) by up to 250 m and the two major outlet glaciers that drain south into Croker Bay (basins 42 and 40) which have experienced advance of between 100 and 400 m along their termini.
Analysis
Variability in the areal changes observed could be due to regional differences in climate forcing, to variations in the sensitivity of glacier mass balance to given changes in climate, or to differences in the response time of different sectors of the ice cap to mass-balance changes. It is also probable, however, that the variations in ice-margin fluctuations have resulted from a complex interaction of these factors.
PREVAILING CLIMATE PATTERN
Differential changes at the ice margin may reflect growth related to the influence of a prevailing climate pattern on the ice-cap geometry. The predominant moisture supply for the ice cap is derived from the North Open Water (NOW) polynya that is located at the north end of Baffin Bay. Cyclonic weather systems that develop over this region deposit over two times the amount of precipitation over the east sector of the ice cap than the west (Koerner, 1977a). Ice-cap growth in the extreme southeast is also encouraged by the maritime effect associated with the NOW, which suppresses summer melt in this region (Koerner, pers. comm., 2002). Koerner (1977a) has suggested that this pattern has likely persisted since at least the Climatic Optimum period (6000–5000 yr BP) to produce the asymmetric geometry of the ice cap that exists today. The eastern basins are typically larger and descend to lower elevations than basins in the west. This results in the east side of the ice cap being subjected to greater amounts of surface ablation and mass loss due to iceberg calving than basins on the west side. Significant mass losses that have occurred along the east side of the ice cap (in contrast to growth along the west) therefore suggest that regions of the ice cap that have experienced preferential growth in the past as a result of this prevailing climate pattern appear to be the most susceptible to increases in average temperature.
BASIN HYPSOMETRY
Analysis of the hypsometry of the major basins was performed to investigate the influence of basin topography on ice-cap response to recent climate changes. Two parameters were used to characterize the basin topography: a Balance Ratio (BR; Furbish and Andrews, 1984) and a sensitivity index derived in this study which provides a measure of the relative change in glacier mass balance resulting from a prescribed change in ELA.
The balance ratio (BR) was calculated from the equation;
where Zm = maximum elevation of the glacier, Zt = terminus elevation of the glacier, and ELA = equilibrium-line altitude. For a rectangular glacier, the change in the altitude of the glacier terminus that results from a prescribed change in the ELA is an inverse linear function of the BR. For glaciers with other shapes, the relationship is still inverse, but nonlinear (Furbish and Andrews, 1984: 202–203). Thus if the ice cap as a whole were subjected to a spatially uniform change in ELA, the change in terminus position resulting from that ELA change would primarily be a function of the BR and planimetric shape of the basin.
BRs were calculated for the major drainage basins of the Devon Ice Cap using ELA estimates derived from the 1960–1966 mass-balance observations by Koerner (1970). The ice cap was divided into quadrants, with boundaries following the watershed divides and appropriate ELA values were applied to each region. ELAs for the northwest and southeast zones were identified directly from Koerner's measurements (950 m a.s.l. and 800 m a.s.l., respectively) and ELAs within the northeast and southwest zones were taken as the average of the two adjacent zones (875 m a.s.l.). Basin hypsometry was derived from the CDED data and analysis was limited to basins larger than 170 km2, which comprise more than 85% of the total ice-cap area.
Although the ELA used to calculate the balance ratio index should be that for an ice cap in steady state, the mass balance associated with the ELAs chosen above is slightly negative. On the basis of analysis of the mass balance/ELA relationship, Koerner (1970), however, suggests that a steady state ELA for the northwest sector of the ice cap is 920 ± 80 m. As the ELA assumed for the northwest region of 950 m falls in this range, we are therefore confident that our method of estimating the ELA for the rest of the ice cap is adequate for proper balance ratio calculations.
The calculated BRs display a strong logarithmic relationship (R2 = 0.78) with percent area change of the major drainage basins (Fig. 6). Most basins with BR values below 1.0 have experienced a significant reduction in surface area. Maximum elevations of these basins are collectively among the lowest on the ice cap (Table 2) resulting in a relatively small altitudinal range above the ELA. In particular, basins 39 and 48 reach maximum elevations of only 998 and 1169 m a.s.l., respectively, resulting in extremely low BR values for these basins. Basins with BRs greater than 1.0 have experienced much less of an area decrease, and in some cases, these basins have increased in size. Those that have experienced growth (42, 40, 25, and 17) are concentrated mainly along the western margin. Although these basins only reach 1500–1600 m a.s.l., their average terminus elevation is 400 m a.s.l., which results in a small elevational range below the ELA. Also, because the termini of these basins are relatively high, they experience cooler temperatures that inhibit glacier retreat. Basins that have experienced limited shrinkage (< −1%) and have BRs greater than 1.0 (2, 44, 41, 43 and 45) extend up to elevations equal to or close to those of the ice-cap summit. These basins generally have larger accumulation area ratios (AARs) than glaciers that have shrunk significantly, and their mass balance was presumably less strongly affected by the post-LIA rise in ELA.
Excluding the southwest arm, basin shape appears to have had some influence on the amount of area loss at the ice-cap margin. Basin shapes were grouped into three main categories according to the Furbish and Andrews (1984) classification scheme: Type “A”—rectangular shaped basins, type “C”—basins that taper up-glacier, and type “D”—basins with wide mid sections that narrow in both up and down glacier directions. For glaciers with BRs less than 1.0, type “D” shaped basins have decreased in size by an average of 3.4% whereas those in type “A” basins have only shrunk by an average of 1.7%. Areal change of glaciers in basins with balance ratios greater than 1.0, show no correlation with shape, however.
Basin hypsometry was analyzed further to determine the sensitivity of major drainage basins to climate warming. Basin sensitivity (BS) values were estimated as the fractional change in AAR resulting from a 100-m rise in the ELA:
where
AAR (accumulation area ratio) represents the proportion of the basin area within the accumulation zone and ΔAAR+100m is the change in accumulation area ratio resulting from a 100-m rise in ELA as documented in the CDED DEM for 1959/1960 (year of the aerial photography used to derive the CDED DEM).
Overall, the calculated BS values display a strong linear relationship (R2 = 0.79) with percent area change of the major drainage basins (Fig. 7). This relationship is controlled mainly by BS values of the basins situated at lower elevations (i.e. basins 38, 52, 48, and 39). These basins have relatively high BS values (>0.2) and have experienced significant areal change over the past 40 yr. The trend is less obvious for basins with BS values <0.2.
The spatial distribution of BS values highlights the relative susceptibility of various regions of the ice cap to climate warming (Fig. 8; Table 2). Basins with BS values greater than 0.23 (“VH”; basins 38, 39, 48, 52) have small AARs and have experienced the greatest average retreat of all major basins since 1960. These basins are restricted to lower elevations (i.e. have low maximum and terminus elevations; see Table 2) where shrinkage appears to be a direct function of the relative loss of accumulation area due to raising the ELA. Drainage basins classified as “highly sensitive” (“H”; 31, 29, 43, 45, 17, and 25) include those that occupy the western margin and have experienced recent advance. These basins are predominantly type “C” in shape implying that the accumulation area ratio increases rapidly as the ELA is lowered. Recent advance throughout this region indicates that these basins may not yet have reacted fully to recent climate warming. Medium sensitivity basins (“M”; 46, 15, 44, and 41) are concentrated on the eastern side of the ice cap. They are predominantly type “A” in planimetric shape (accumulation area ratio decrease is directly proportional to a rise in ELA) and reach near-summit elevations with AARs averaging 0.71. Finally, low sensitivity drainage basins (“L”; 12, 42, 2, and 40) have properties very similar to the “M” basins except that they drain towards the north and south coasts and have an average AAR of 0.81.
The hypsometric curves plotted in Figure 9 indicate the distribution of surface area within a basin in relation to elevation above or below the ELA. Basins in which the ELA intersects the steepest section of the curve (basins 38 and 52) are most sensitive to climate change because moving the ELA from the current location results in changes in accumulation area that are large relative to the total basin area. In some cases however, the ELA intersects the hypsometric curve close to the upper limit of the cumulative basin area (basins 39 and 48). The accumulation areas within these basins have almost entirely disappeared resulting in consistently high ELAs regardless of fluctuations of the ELA. The ELAs of highly sensitive basins (“H”) intersect steep sections of the hypsometric curves. The sensitivity of these basins is reduced however by the fact that the cumulative area below the ELA of these basins is generally less than 0.4. Medium-sensitivity basins (“M”) have between 0.2 and 0.3 of their cumulative area below the ELA whereas low-sensitivity basins consistently have less than 0.2 of their area below the ELA. Both of these basin types intersect relatively low sloping sections near the inflection point along the hypsometric curves indicating that the accumulation area of these basins will decrease at a much faster rate as the ELA is raised than it will increase as the ELA is lowered.
RESPONSE TIME
Approximate response times (RT) were calculated for each of the major drainage basins using the volume time scale formula in Johannesson et al. (1989):
where H is the maximum thickness of ice within the basin and bt is the mass balance at the glacier terminus. Ice thickness data were obtained from Dowdeswell's 2000 RES data set (Dowdeswell et al., in press) and terminus ablation rates were interpolated from Koerner's, 1961–1998 mass-balance observations as compiled by Dyurgerov (2002).
A moderate inverse relationship (R2 = 0.65) exists between RT and percent area change of the major basins across the ice cap (Fig. 10). This relationship indicates that maximum areal changes have occurred for basins with the shortest calculated response times, and that the magnitude of change decreases logarithmically (and in some cases the basins grow rather than shrink) as the calculated response times increase.
Three distinct groupings emerged from these RT calculations (Table 3). RTs in the first group of basins range from ca. 125 yr to ca. 380 yr. This group includes ∼70% of all major basins. Although all basins in this group terminate at or near sea level where ablation rates are highest, basins with the shortest RTs in this group (29, 38, 39, 46, 48, and 52) are restricted to the peripheral regions (i.e., their source regions do not extend to the ice-cap summit region) where the ice is thinnest. Except for basin 39 (southwest arm), all basins in this category are located along the east or south margins of the ice cap. Other basins within this group (2, 40, 42, 12, 41, 15, and 44) extend farther inland towards the ice-cap summit, lengthening the RT due to increased ice thickness. Basins exhibiting intermediate (ca. 700 yr) and long (ca. 950 yr) RTs terminate at higher elevations (500–600 m a.s.l.) where terminus ablation rates are relatively small. These basins are concentrated along the western margin and they also extend inland towards the interior summit region where ice thickness is greatest.
The RTs calculated here are consistent with characteristic RTs suggested by Bahr et al. (1998). Glaciers in larger basins with steeper mass-balance gradients (east side) have significantly shorter RTs than glaciers in smaller basins with shallower mass-balance gradients (west margin). This difference suggests that the east and west sectors of the ice cap may be responding to distinctly separate periods of climate forcing.
MULTIPLE REGRESSION ANALYSIS OF AREAL CHANGE
To determine the influence of BS, BR, and RT on the areal changes (AC) measured, multiple regression analysis was next performed resulting in the equation:
BR was highly correlated with RT (r2 = 0.86) and was therefore excluded from this analysis. P-values of 0.00001, and 0.002, were obtained for the coefficients of BS and RT, respectively. The coefficient of multiple determination (R2 = 0.854) indicated that 85.4% of the total variance of AC is explained by BS and RT.
VOLUME CHANGE ESTIMATES
Volume change was calculated for the ice cap using two independent methods. The first method was based on calculating the difference between the 1960 ice-cap volume, derived from the volume-area scaling (VA) technique (Bahr et al., 1997), and the 2000 ice-cap volume derived from Dowdeswell et al.'s (in press) RES data. Using this method, the 1960s volume was calculated for the southwest arm and the main ice cap separately, and summed for a total ice-cap volume estimate. The second technique was the method presented by Hooke (1998: 219), which calculates the change in volume of a glacier as the product of the maximum thickness and the length change at the terminus (MT). In order to account for the influence of glacier width, however, the total area of change measured in each basin was multiplied by the maximum thickness to estimate the resulting change in volume. Using this method, volume change was calculated for each major drainage basin separately and summed to estimate the total ice-cap volume change.
The 1960 ice-cap volume (VOL1960) was calculated from the VA method using the formula:
where, A1960 is the area of the ice cap in 1960 and 1.25 is the appropriate exponent for small ice caps (Bahr et al., 1997). The scaling constant (SC) is derived from the formula: where A1999 is the area in 1999 and VOL2000 is the volume in 2000 taken from Dowdeswell et al.'s (in press) 2000 RES data. VOL2000 measurements of 3854 km3 and 130 km3 produced scaling constants of 0.031525 and 0.009968 for the main ice cap and the southwest arm, respectively. Applying the scaling constants to the 1960 measurements resulted in a total volume estimate of 4181 km3 (146 and 4034 km3 for the southwest arm and main ice cap respectively). Error was estimated for the VA method at ±12 km3 based on the 1960s area measurement error of ±30 km2. Error for the MT method of ±9 km3 was calculated as a product of the average area measurement error between 1960 and 1999 (±35 km2) and the average ice-cap thickness (284 m).
Differences between the 1960 and 2000 volume estimates derived from the VA method and from Dowdeswell et al.'s (in press) RES measurements respectively indicate a loss of 71 km3 for the entire ice cap and 16 km3 for the southwest arm over the past 40 yr. Volume changes derived from the MT method were similarly estimated as a loss of 62 km3 for the ice cap in total and 18 km3 over the southwest arm. Total volume loss of the Devon Ice Cap between 1960 and 1999 was therefore estimated as 67 ± 12 km3. These estimates of volume loss imply a thinning of ∼3 m averaged over the main part of the ice cap and ∼8 m over the southwest arm.
Our computed volume loss of 67 km3 equates to a 0.17 mm contribution to global sea level between 1960 and 1999. When extrapolated to all ice cover in the Canadian Arctic, this would represent a total contribution of ∼1.7 mm. This is approximately one-sixth the current estimated contribution from all Alaskan glaciers (Arendt, et al., 2002) which cover ∼40% less surface area than ice in the Canadian Arctic.
Conclusions
Remotely sensed imagery acquired in 1959/1960 and 1999/2000 reveals a net decrease in surface area of the Devon Ice Cap of 332 ± 40 km2 with an associated volume loss of 67 ± 12 km3. The dominant trends in geometric change currently being experienced by the Devon Ice Cap are retreat (thinning) of the eastern margin, growth of the western basins, and rapid shrinkage of the southwest arm. Major tidewater glaciers along the east coast have retreated up to 1.3 km over the past 40 yr. These glaciers drain from basins that have experienced increased bedrock exposure indicating that the ice surface within these basins has lowered. Advance of the northwest margin suggests either that conditions favourable to ice sheet growth are prevalent in this region, or that this sector of the ice cap is still responding to the cooler conditions of the Little Ice Age. The southwest arm of the Devon Ice Cap lies almost entirely below the regional ELA (Koerner, 1970). Thus, rapid retreat of the ice margin in this region reflects the loss of the accumulation area for this sector as a result of ELA rise driven by post-LIA warming.
Statistical analysis has shown that 85.4% of the variance in basin area change is explained by variations in the balance ratio, basin sensitivity, and response times of different drainage basins. The logarithmic relationship between BR and AC suggests that the ice cap will shrink at a much faster rate in response to a given rise in the ELA than it will grow in response to an equivalent lowering of the ELA. Similarly, cumulative area curves for most basins indicate that the rate of accumulation area loss as the ELA is raised is greater than the rate of accumulation area gain when the ELA is lowered. Finally, response times of 700 to 1000 yr throughout the western region indicate that the advance of these basins is likely a delayed response to cooler LIA climatic conditions. In contrast, shorter response times (100–400 yr) for the remaining basins suggest that the observed retreat over these regions is likely a response to post-LIA warming.
Acknowledgments
Financial support for this research was provided by an NSERC Discovery grant to Dr. Martin Sharp. Funding for David Burgess was provided through an NSERC IPS-2 scholarship sponsored by Land Data Technologies Inc. of Edmonton, Canada. Landsat 7 ETM+ imagery and Canadian Digital Elevation Data (CDED) were purchased through support from CRYSYS (Cryospheric System to Monitor Global Change in Canada), the Meteorological Services of Canada (MSC) and the Canadian Space Agency (CSA—Climate and Cryosphere Initiative). Aerial photography was provided by the National Hydrology Institute and Art Dyke of the Geological Survey of Canada. We thank Prof. J. Dowdeswell for making available radio echo sounding data for the Devon Ice Cap, Dr. John Andrews and an anonymous reviewer for their comments to the manuscript, and Dr. R. M. Koerner for his comments on an earlier draft.
References Cited
FIGURE 1.
Location of the Devon Island Ice Cap in the Canadian Arctic Archipelago. Image of the ice cap. is an orthorectified mosaic of Landsat 7 satellite panchromatic data acquired on 13 July 1999 and 2 August 2000. Contour interval is 200 m
FIGURE 2.
Drainage basins of the Devon Ice Cap. Only the basins labeled with unique identification numbers were used for analysis in this study
FIGURE 3.
1960 ice margins superimposed on 1999/2000 Landsat image mosaics highlight retreat of the major tidewater glaciers over the past 40 yr. Shaded basins have experienced a minimum increase of 2 km2 of exposed bedrock area between 1960 and 1999 within the ice cap interior
FIGURE 6.
Plot of the logarithmic relationship between calculated balance ratios and percent area change of glaciers in the major drainage basins between 1959/1960 and 1999/2000
FIGURE 7.
Plot of the linear relationship between basin sensitivity and percent area change of glaciers in the major drainage basins measured between 1959/1960 and 1999/2000
FIGURE 8.
Relative sensitivity of glaciers in the major drainage basins to increasing the ELA by 100 m
TABLE 1
Summary of error associated with extracting linework from the satellite and aerial photography image data
TABLE 2
Hypsometrically derived index values and associated glacier attributes for all major drainage basins. Shaded rows indicate basins with balance ratios below 1.0
TABLE 3
Response time estimates and associated glacier attributes for each major drainage basin across the ice cap. Shaded rows highlight the intermediate class of response time estimates
